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Abstract 

An overview of Hamiltonian systems with noncanonical Poisson 
structures is given. Examples of bi-Hamiltonian ode's, pde's and lat- 
tice equations are presented. Numerical integrators using generating 
functions, Hamiltonian splitting, symplectic Runge-Kutta methods are 
discussed for Lie-Poisson systems and Hamiltonian systems with a gen- 
eral Poisson structure. Nambu-Poisson systems and the discrete gra- 
dient methods are also presented. 

Keywords: Hamiltonian ode's and pde's, symplectic integrators, 
Lie-Poisson systems, bi-Hamiltonian systems, integrable discretizations, 
Nambu-Hamiltonian systems. 

1 Introduction 

In this review article we present Hamiltonian systems with a non- 
canonical Poisson structure and give an overview of Poisson structure 
preserving numerical integrators. In the last decade, efficient and re- 
liable numerical integrators were constructed for canonical Hamilto- 
nian systems with symplectic structure (see The theoretical 
foundation of stability and convergence of symplectic integrators was 
provided by the backward error analysis j^H- Methods based on the 
preservation of the geometric features of the continuous system by dis- 
cretization are called "geometric integrators", which include also the 
symplectic methods (see for a state of art review of these methods and 
their applications |10i). They have better stability properties, smaller 
local and global errors than the standard integrators based on the 
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local error control. The quadratic invariants of the underlying contin- 
uous system are exactly preserved by symplectic integrators. Compli- 
cated invariants, among them the polynomial invariants of order three 
and greater are in general not preserved (^01) USD- The phase space 
structure of low dimensional Hamiltonian systems are very accurately 
preserved by symplectic methods. But the applicability of symplectic 
integrators for higher dimensional systems arising from semidiscretiza- 
tion of Hamiltonian pde's is not well established 1. 

The noncanonical Hamiltonian systems on Poisson manifolds arise 
in many applications. The best known example is the Euler top equa- 
tion, a three dimensional Hamiltonian system. The Poisson structure 
is more general than the symplectic structure, the dimension of a Pois- 
son manifold doesn't need to be even as in the case of the Euler top 
equation compared to symplectic manifolds which appear only in even 
dimensions. The Poisson systems occur as finite dimensional systems 
in the form of ode's and nonlinear lattice equations or as infinite di- 
mensional systems in the form of pde's. Hamiltonian pde's like the 
Korteweg deVries (KdV) equation, the nonlinear Schrodinger equation 
(NLS) and Euler equations of an incompressible fluid have many ap- 
plications in hydrodynamics, optics and meteorology. The algebraic 
and geometric properties of the Poisson manifolds, like the Poisson 
bracket, Darboux transformations, Casimirs, symplectic foliation and 
integrability, are given in Section 2. Many Hamiltonian systems can 
be written with respect to the same coordinate system in more than 
one Hamiltonian formulation. The so called bi-Hamiltonian systems 
occur as ode's, pde's and in nonlinear lattices. The properties of bi- 
Hamiltonian systems, and finite dimensional examples are given in 
Section 3. The real application of the bi-Hamiltonian formulation is 
the field of Hamiltonian pde's and nonlinear lattice equations, which 
are presented in Section 4. 

Geometric integrators for noncanonical Hamiltonian systems or 
Poisson systems are largely unexamined. There exist a variety of meth- 
ods for the Lie-Poisson systems, Hamiltonian systems with a linear 
Poisson structure. The generating function approach ^31 and splitting 
methods based on the separability of Hamiltonian for canonical Hamil- 
tonian systems were extended to the Lie-Poisson systems |M| ■ Unfor- 
tunately there exits no general method for Hamiltonian systems with 
a general Poisson structure beyond the Lie-Poisson integrators. One 
can show that some of the symplectic methods for canonical Hamilto- 
nian systems may preserve the Poisson structure of some Hamiltonian 
systems. This means that whether an integrator is Poisson or not, de- 
pends on the specific problem. For some nonlinear lattice equations, it 
is possible by using the Darboux coordinates to transform the Hamil- 
tonian system with the Poisson structure to a canonical one and then 
integrate this system by the symplectic methods. But there exists no 
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general algorithm for such a transformation. These aspects of Poisson 
integrators are discussed in Section 5. 

In Section 6 the so called Nambu-Hamilton systems which can be 
considered as an extension of the Poisson structure are presented. The 
structure of these systems is discussed and the discrete gradient meth- 
ods which applicable for them, are described. 

2 Poisson structures and Hamiltonian sys- 
tems 

In order to present the finite dimensional and infinite dimensional 
Hamiltonian systems, we have to look at the properties of the Poisson 
manifolds and functions defined on these manifolds. In our presenta- 
tion we follow mainly Chapters 6 & 7. The key property here is 
the so called Poisson bracket, which assigns to each pair of functions F 
and G in the space of smooth functions ^{M) on a smooth manifold 
Ai a third function, denoted by {F, G}. 

Definition 1 The Poisson bracket is a bilinear operation satisfying 
the following conditions: 

(i) skew- symmetry 

{F,G} = -{G,n, 

(ii) the Leibniz rule 

{F,G-P} = {F,G]-P^G-{F,P], 
(ii) the Jacobi identity 

{{F, G}, P} + {{P, F}, G} + {{G, P}, F} - 0, 

where ■ denotes the ordinary multiplication of functions. 

A manifold Ai with a Poisson bracket is called a Poisson manifold and 
the bracket defines a Poisson structure on A4. The Poisson manifold 
is more general then a symplectic manifold; it doesn't need to be even 
dimensional. The canonical Poisson bracket on a even dimensional 
Euclidean space A4 = M™, ni — 2n with the canonical coordinates 
{p, q) = (pi, . . . ,p„, gi, . . . , (7„) is given for two smooth functions F(p, q) 
and G{p, q) by 

with the bracket identities 

{Pi,P]} = ^, {9i,9i} = 0, {qi,pj} = 6], i,j = l,...,n 
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where Sj denotes the Kronecker symboL 

The general local coordinate picture of a to dimensional Poisson 
manifold is given by the Poisson bracket on local coordinates x = 
{xi, . . . , x„i). The Poisson bracket for two smooth functions F{x) and 
G{x) takes the form 

(f.Gl.fffe..,}^^ (2, 

i=l j=l -I 

A Hamiltonian vector field Xh on AA with the coefficients functions 
S,i{x) depending on H is given by 

Z— 1 

and satisfies 

Xh{F)^{F,H} = -{H,F}. 

The governing equations of the flow associated with the vector field 
are called Hamiltonian equations 

x^{x,H}. (3) 



Definition 2 A general Hamiltonian system is given then by the triple 
M., the manifold, {•, the Poisson structure and H{x) a real function 
on M. The Poisson structure gives the coordinate functions Jij = 
{xi,Xj} which are called structure functions. They satisfy the Leibniz 
rule and have the following properties 

i) skew- symmetry: 

Jij{x) ^ -Jji{x), i,j = l,...,m 

ii) Jacobi identity: 

m 

^{JiidiJjk + JkidiJij + JjidiJki} = 0, i,j = l,...,m (4) 

i=l 

for all X € M with di = d/dxi. 

The conditions for the Jacobi identity Q form a system of nonlinear 
partial differential equations for the structure functions Jij{x). The 
structure functions can be assembled into a skew-symmetric structure 
matrix J{x) G R™^™. In this case the Poisson bracket (0) for the 
functions F and G takes the form 

{F, G] = VF{xfj{x)VG{x) 
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and the Hamiltonian equation Q can be rewritten as 

X = J{x)VH{x). 

The most important examples of Poisson structures are associated with 
m-dimensional Lie algebras. The associated Lie-Poisson bracket be- 
tween two functions F,G with the structure constants c^^ ihjjk = 
1, . . . , TO of the Lie algebra, is given by 

j,k—l •' 

with the linear structure functions Jij{x) = X^fcLi cltj^k- 



Example 1: Euler top equation 

The most important example of a Lie-Poisson system 

dxi I, dH 

is the Euler top equation describing the motion of a rigid body around 
a fixed point: 

I2 — h ■ Iz — h ■ h — 12 /r\ 

Xl = —— 2:22:3, X2 = 2:3X1, X3 = — — 0:1X2 (5) 

hJ^i ^3^1 ^1^2 

where /i,/2,/3 denote the moments of inertia. The Poisson structure 
is given by 

{xi,X2} = -X3, {X2,X3} = -2:1, {X3,xi} = -X2 (6) 

It can be written in the Hamiltonian form (jSJ with the energy integral 
H as the Hamiltonian 

and with the linear structure matrix 

(0 — X3 X2 
X3 — Xl 
— X2 Xl 

cIt 

— = J{x)VH{x) = X X VH{x) 



here x denotes the cross product. 



5 



Theorem 1 A Poisson manifold of dimension m is symplectic if the 
rank of the Poisson structure is maximal everywhere. Because the rank 
of a Poisson manifold at any point is always an even integer, a sym- 
plectic manifold is necessarily even dimensional. In other words if the 
structure matrix J{x) is invertible, then according to Darboux's the- 
orem the Poisson structure is isomorphic to the standard symplectic 
structure in R^". In this case the structure matrix J{x) determines a 
symplectic structure on A4 C M™, m = 2n if and only if its inverse 
K{x) = J{x)~^ satisfies the following conditions (14 If, Proposition 
6.15) 

i) skew- symmetry: 

Kij{x) ^ -Kji{x), i,j = l,...,m 

ii) Jacobi identity: 

dkKij{x) + djKki{x) + diKjk{x) =0, i,j,k^l,...,m 

The nonlinear partial differential equations for the structure functions 
in the Jacobi identity ^ simplify to linear partial differential equations 
involving Kij (x) . 

Definition 3 gT^ A Poiss on structure is called degenerate, if the struc- 
ture matrix J{x)is non- invertible. It follows from the skew-symmetry 
that the rank of the Poisson matrix at any point is always an even 
integer, which implies that odd dimensional systems have a degenerate 
Poisson structure. If the manifold is of dimension m = 2n-\- I and the 
rank of J{x) is 2n everywhere, then there exist distinguished functions 
Fi{x) = Ci,i = 1, . . . ,1, called Casimirs, whose Poisson bracket vanish 
with any function or variable Xi , i.e. {Fi{x), Xi} = Q,i = \, . . . ,1. 

If the Poisson structure is of constant rank, the symplectic foUation 
simpHfies. One can introduce local coordinates, which bring the foli- 
ation to the canonical form, which results from the Darboux theorem 
(|41j. Theorem 6.22). At each point a; G of an m-dimensional man- 
ifold of constant rank 2n < m everywhere, the Darboux theorem says 
then, there exist local coordinates (gi, . . . , g„,pi, . . . ,p„, zi, . . . , z;) in 
terms of which the Poisson bracket takes the canonical form where 
Pi,qi denote the canonical variables. There exist locally I indepen- 
dent Casimirs Fi{x), ■ ■ ■ , Fi{x); their common level surfaces Fi{x) = 
Ci, - ■ ■ ,Fi{x) = ci are called symplectic leaves of the Poisson bracket. 
Using the Darboux transformation, the Poisson structure can be trans- 
formed on the symplectic leaves to the canonical one. 

If the structure matrix J{x) is constant and nonsingular in some 
coordinate systems, then a linear change of coordinates reduces the 
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Poisson structure to the canonical one JTj . The reduction of the Pois- 
son structure using Casimirs on symplectic leaves produce usually more 
complicated systems. For example for the Euler top equations, the to- 
tal angular momentum 

1^1(^1 + 4 + 4). (8) 

is a Casimir (commutes with Xj, i = 1, . . . , 3) and the symplectic leaves 
are the spheres described by ©■ Using the transformations 

the Euler top equation can be reduced to a canonical Hamiltonian 
system (0, pp. 7) 

H{q,p,c) = {I -p^) f — cos2(7+ — sin^gj + — . 

The geometric characterization of the local and global structure of 
Poisson manifolds is given by the splitting theorem of Weinstein |221 : 

Theorem 2 On a Poisson manifold Ai, any point x ^ JV[ has a co- 
ordinate neighborhood with coordinates qi, . . . ,qk,pi, ■ ■ ■ ,Pk,yi, ■ ■ ■ ,yi 
centered at x, such that 

This theorem gives a decomposition of the neighborhood of x as the 
product of two Poisson manifolds: one with rank 2k and the other with 
rank at x. The local structure of a Poisson manifold is described 
by the symplectic leaves and their transverse Poisson structures (the 
second part in\^. Globally, a Poisson manifold is obtained by gluing 
together the symplectic manifolds. There are three special cases of the 
splitting theorem i5^ : 

• If the rank is locally constant, then fij = 0. The splitting theorem 
reduces to the Lie 's theorem, i. e. each point of A4 is contained 
in a local coordinate system with respect to which the structure 
functions Jij are constant 

• At the origin of A4 one has only yi 's, the first term in ^ does 
not appear. 

• A symplectic manifold is a Poisson manifold where rank Y[ — 
dim M. everywhere. The splitting theorem gives canonical coor- 
dinates gi, • • • • ■ • ,Pk- 
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There exists a close relation between the symmetries and conservation 
laws of Hamiltonian systems. Noether's theorem provides a connection 
between the symmetries and conservation laws or first integrals and can 
be applied to Hamiltonian systems. The first integrals are character- 
ized by the vanishing of the Poisson bracket: a function X is a first 
integral of the autonomous Hamiltonian system (PJ, if {X, H} — for 
all a; G A^. A Hamiltonian system with a time- independent Hamil- 
tonian has H{x) as the first integral. The Casimirs Fi{x) determined 
by the Poisson structure J{x) are first integrals too, but they arise 
from the degeneracies of the Poisson bracket and are not generated by 
the symmetry properties of the Hamiltonian system. If the Poisson 
manifold is symplectic, only constant functions are Casimirs. 

Hamiltonian systems are usually non-integrable, there are in gen- 
eral no additional integrals of motion. But there is a class of Hamilto- 
nians which are called completely integrable. 

Definition 4 141^ A Hamiltonian system is called completely inte- 
grable if it has n integrals of motion Ik{q,p),k = 1, . . . ,n 

• which do not depend explicitly on time, 

• are functionally independent, i.e. 

• are in involution 

{Ij, Ik} = for each pair {j,k), i,k — 1, . . . , n. 

The Liouville theorem then says that any integrable system can be 
integrated by quadratures. The geometric content of complete inte- 
grability is described by the Liouville- Arnold theorem: 

Theorem 3 \44\ If ihe surface of the constant value of the 
Hamiltonian H{q,p) ~ E is compact connected manifold, then 

• the manifold 

M = {{q,p) : Ikiq,p) = Ck = const, k = l,...,n} 

diffeomorphic to the torus T", 

• in the canonically conjugate action and angle coordinates (Ik, ^k), 
Hamilton equations have a simple form: 

ipi = uJk{I), ik = 0, k = l,...,n, 
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• the action- angle variables have a simple dependence on time 
hit) const, ipk(t) = (^[°' + LOkt, 



the variable ipk determines the position of a point on the torus 
and vary linearly with time; the motion is quasiperiodic on the 
torus T", 

• the Hamilton equations can be solved by quadratures. 

The Euler top equation is completely integrable. Other examples of 
completely integrable rigid body equations are the Lagrange top, and 
the Kowalevskaya top equations 0] . 

If a Hamiltonian system has 2n — 1 independent integrals at most, 
it is called completely degenerate. Such systems occur in very rare 
cases for Hamiltonian with a potential function V{q) = a/q for the 
Kepler problem, V{q) = ^q^ for the harmonic oscillator and for the 
Calogero-Moser system |48| . 

3 Bi-Hamiltonian systems 

Definition 5 A system of differential equations is called bi-Hamiltonian 
if it can be written in the Hamiltonian form in two distinct ways ■' 



with the structure matrices Jo{x) and Ji{x) determining the Poisson 
bracket {F,G}k = Vi^"^ Jfc(x)VG, fc = 0, 1. The structure defined by 
Jo{x) and Ji{x) is called a Hamiltonian pair. 

A Hamiltonian pair is compatible, if the linear combination of Ji{x) — 
XJo{x) determines also a Poisson bracket for any real A. The compat- 
ibility condition is not trivial, because the differential equations must 
be satisfied for the Jacobi identity I^J. In the symplectic case (non- 
degenerate Jk{x), k — 1,2 and n is even), the Jacobi conditions can be 
replaced by the symplectic two-forms 



which are closed, i.e. dQk — 0. For a given pair, the corresponding 
Hamiltonian system can be found by solving the linear system of partial 
differential equations Vi?i = J^^ {x)Jo{x)S/Hq ^U]. A Hamiltonian 
pair is called non-degenerate at the point x if the skew-symmetric ma- 
trix pencil J\ (x) = Ji (x) — A Jo (a;) is nonsingular for at least one finite 
A. According to Magri bi-Hamiltonian systems with a compatible 
Hamiltonian pair are completely integrable. 



— = Jo(x)Vi/i = Ji{x)VH„ 



(10) 
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Theorem 4 For any bi-Hamiltonian system with a non- degenerate, 
compatible bi-Hamiltonian structure, there exists a hierarchy of Hamil- 
tonian functions Hq, Hi, H2, ■ ■ ., all in involution w.r.t. both Poisson 
brackets {Hi, Hj}k = 0, k — 1,2, i, j — 0, . . . and generating mutually 
commuting bi-Hamiltonian flows \4^ ), they are completely 
integrable. The sequence of first integrals can be constructed for bi- 
Hamiltonian systems by the Lenard recursion scheme (J^ [3 \^): 

{■,Hq}o = Hq is a Casimir with respect to { , }o 
{■,Hi}o = {■,Ho}i, 
{■,H2}o = {■,Hi}i, 

{■,Hk+i}o = {■,Hk}i 
which generate a sequence of bi-Hamiltonian vector fields 

Xk = {x,Hk}Q, 

= {x,Hk^i}i, A: = 1,2,... 

In case of infinite dimensional systems (pde'e), one has an infinite 
number of first integrals. For finite dimensional systems, the Lenard 
recursion is either finite or infinite, in the second case only a finite 
number of integrals are linearly independent 9 . 

Example 2: Euler top equation as bi-Hamiltonian system 

The Euler top equation Q is also bi-Hamiltonian with respect to the 
second Poisson bracket 

{xi,X2}^^, {x2,X3}^^, {2:3,2:1} = ^ (11) 

J3 J2 

with the total angular momentum JSJ as the second Hamiltonian. The 
corresponding structure matrix is given by 



Ji{x) 



Examples of finite dimensional Hamiltonian systems with Poisson 
and bi-Hamiltonian structures are given in |25 |. |44|. 

Example 3: Lotka-Volterra equations 

Three dimensional Lotka-Volterra system of competing species has a 
bi-Hamiltonian structure 

-2;i2;2 2:i2;3 \ / X1X2X3 -X1X2X3 

X1X2 — 2:22:3 V-ffi = — 2:12:22:3 2:12:22:3 

—2:12:3 2:2X3 y y 2:12:22:3 —2:12:22:3 

(12) 
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with the Hamiltonians 

Hi = xi + X2 + X3 and H2 = log(a;i) + log(.X2) + log(a;3). (13) 
Example 4: Lorenz system 

xi = a{x2 - xi) 

X2 = {pXl - XiXs, - X2) (14) 
2^3 = {-I3x^+ X1X2) 

has a chaotic behavior for most values of the parameters. It is Hamil- 
tonian with the following parameter values and the Hamiltonian 

1 _T? -I- .T? 

p = 0, cr = -, /? = 1, H 



2"- ' {xl-X:,Y- 
with respect to the Poisson structure 

{xi,X2} ^ xfx^+xl, {a;2,a;3} = -2(a;2+X3)a;i, {x3,xi} ^ {xl-X3)x2. 

The system has two time-dependent conserved quantities. Using a 
transformation of the dynamical variables and time, the equation H14|l 
takes a bi-Hamiltonian form [221 

, 1 

"i = 2^2, U2 = -U1U3, U3 = W1U2. 

with the Hamiltonians Hi = — uf, H2 = ^2 + ^3 ^^'^ with the 
corresponding Poisson brackets 

{ui,-U2}o = ^, {u2,U3}n = -^Ui, {ui,U3}o=0, 
{wi,U2}l = ^U3, {w2,U3}l = 0, {ui,U3}^^U2 

Other three dimensional systems admitting bi-Hamiltonian struc- 
ture like the Halphen system, May Leonard equations, Maxwell-Bloch 
equations, and the Kermack-McKendrick model for epidemics are given 
in [23 El- 
There are few examples of bi-Hamiltonian systems in dimensions 
higher than three W- We give here as an example the Calogero-Moser 
system with two particles. 

Example 5: Calogero-Moser system 

The motion of repelling particles on a line is described by the Hamil- 
tonian 



1 1 
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The four dimensional bi-Haniiltonian Calogero-Moser system for two 
particles is given in |32j with the Hamiltonians Hq and Hi. 



2 [qi - q2r 
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P2 



29^2 

49?2+P?2 
_ ^9l2Pl2 
49^2 +P?2 

29?2 



where gi2 



and pi2 =Pi -P2- 



4 Bi-Hamiltonian pde's and lattice equa- 
tions 

Many nonlinear partial differential equations like the Korteweg de Vries 
(KdV) equation, the Eulcr equations of fluid dynamics and the non- 
linear Schrodinger equation can be written in the Hamiltonian form. 
The Hamiltonian formulation for the ode's can be extended to pde's by 
replacing the Hamiltonian H{x) by a Hamiltonian functional 7i[u], the 
gradient operation Vi? by the variational derivative S7i and the skew- 
symmetric matrix J(x) by a skew-adjoint differential operator 0(2?). 
The resulting Hamiltonian system has the form 

^ - e(P) X SH[u] 

The right hand side of the Hamiltonian evolution equations 

ut = K[u] ^ K{x,u'^'^^) 

contains the space variable vector x and m*-"^ represents all derivatives 
of the function u{x, t) with respect to x of order at most n. The Poisson 
bracket between two functional J- and Q is defined by 

{:F,g} ^ J ST-e{v)gdx 

which satisfies the skew-symmetry condition and the Jacobi identity. 
An interesting feature of the Hamiltonian pde's is that many of them 
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have a bi-Hamiltonian structure. As a result of this, there exists a 
hierarchy of mutuaUy commuting Hamiltonian flows, and conserved 
quantities, generated by the recursion operator based on two compat- 
ible Poisson operators. In the last two decades the bi-Hamiltonian 
structures of many pde's were discovered, like the KdV equation, the 
nonlinear Schrodinger equations, the Euler equation, and the Boussi- 
nesq equation 01] . 

The Korteweg de Vries (KdV) equation has the bi-Hamiltonian 
form 

ut = Uxxx + uux = Dx5 J (^^^^l + ^"^^ = QoSHi{u) 

2 1 \ /" 1 

Dl + -uDx + -UxU / -u^dx = BiSHoiu) 

The nonlinear Schrodinger (NLS) equation with the complex vari- 
able u 

iut + Uxx + 2|u|^u = 0, 
has the following bi-Hamiltonian structure IH] 

Ut \ I iuxx + 2m|?ip 

— m*^ — 2m*|up 



i 



^5 j^{\ux\^ -\ut)dx 



2uD-^u -D-2uD-^u* \ , f I ., ^ 
- 2u*D-^u 2u*D-^u* ) ^J^ 2*^"^" " ^^"^ 

with u* denoting the complex conjugate of u. Both systems, the KdV 
and NLS equations, possess due to their bi-Hamiltonian structure, in- 
finitely many first integrals, i.e. they are completely integrable. The 
integrals can be generated using the Lenard scheme. If the differ- 
ential operators Oq and Qi are compatible, the infinite hierarchy of 
bi-Hamiltonian systems is defined recursively by 

Ut ^ OqSH^+i = BiSHi, i = -l,... 

with the integro-differential recursion operator R = QiO^^ 



Nonlinear lattice equations or discrete nonlinear systems play an 
important role in modelling many physical phenomena. Examples of 
them are the Fermi-Pasta-Ulam lattice, the Toda and Volterra lattices, 
the Calogero-Moser system. They also arise by semi-discretization of 
some integrable nonlinear pde's like the Ablowitz-Ladik lattice which 
represents an integrable discretization of the nonlinear Schrodinger 
equation Most of the nonlinear lattices are completely integrable. 
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The conserved quantities, the multi-Haniiltonian structure, the recur- 
sion operator, Lax form, and master symmetries of the lattice equations 
were investigated in recent years by many authors (see for example 

nsiinni). 

Example 6: Toda lattice 

The motion of n interacting particles with an exponential interaction 
potential Hamiltonian 

n 

are described by the Toda lattice equations. Here the g^s denote the 
positions of the z-th particle and the p^s their momentum. Introducing 
the new variables = 6*^'^'"'^'+'^^ hi = pt after Flaschka transforma- 
tion, one obtains BHl[T?)j 

di ^ ai{bi+i - b,), b, ^ ai - a,^i, l<i<n (15) 

with the open-end (ao = a„ = 0) or periodic (ao — a„, 6„+i = bi) 
boundary conditions. It has a tri-Hamiltonian structure with the fol- 
lowing Poisson brackets and corresponding Hamiltonians j55| 

• linear Poisson bracket (Lie-Poisson bracket): 



{bi,a.i}i = ai, {ai,bi+i}i = ai 

n N 



2 

1=1 i=l 



• quadratic Poisson bracket: 

{h, ai}2 = aibi, {o^, 6i+i}2 = a^b^+i 

{ttj, ai+i}2 = Oi+iOi, 6i+i}2 = flj 

n 
i=l 

• cubic Poisson bracket: 

{bi,ai}3 = a^{bf +Ut), {ui^bi+ijs = atibf+i + a.^) 

{ai,ai+i}3 = 2aiai+ibi+i, {bi,bi+i}3 = a,,{bi + bi+i) 



1 " 

fi = ^ E ^°S(a< 



2 

i=l 
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With respect to the hnear and quadratic Poisson brackets, the Toda 
lattice forms a bi-Hamiltonian system. Because both brackets are de- 
generate, it is not possible to find the recursion operator and apply the 
Lenard scheme to find the hierarchy of the integrals, which was done 
using master symmetries A different bi-Hamiltonian formulation 
of the Toda can be found in ^1 

Another interesting nonlinear lattice equation is the Volterra lattice 
[SS] which is a model for vibrations of the particles on lattices (Liou- 
ville model on the lattice) and describes the population evolution in a 
hierarchical system of competing species. It can be also considered as 
an integrable discretization of Korteweg de Vries equations. 

Example 7: Volterra lattice 

yi = yi{yi+i~yi-i), i = i,---,n (16) 

is bi-Hamiltonian 

y = Jo(y)Vi/i = Ji(y)Vifo 
with respect to the quadratic and cubic Poisson brackets 

{yi,yi+i}a ^yiVi+i, (17) 
{yz,yt+i}i = ytVt+iiyz + 2/,:+i), {yt,yt+2}i = ytyt+m+2- (18) 

and the corresponding Hamiltonians are 

n 1 ^ 

i=l 1=1 

The Toda and Volterra lattices are closely related through the nonlin- 
ear transformation of variables ^Hl 

Oi = 2/2i2/2j-i, h = y2i-i - y2i-2- 

5 Lie-Poisson and Poisson integrators 

It is natural to ask if the symplectic integrators applied to canonical 
Hamiltonian systems on symplectic spaces can be extended to produce 
integrators that preserve the linear Poisson structure (Lie-Poisson in- 
tegrators) or the more general Poisson structures. 

5.1 Integrators based on generating functions 

The generating function approach in the framework of symplectic inte- 
grators was applied by several authors to Lie-Poisson systems. Ge and 
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Marsden constructed Lie-Poisson integrators in |24j by finding approxi- 
mate solutions to the Hamilton-Jacobi Lie-Poisson equations using the 
generating functions developed in (13) for canonical Hamiltonian sys- 
tems which is symplectic and energy preserving. The time- dependent 
Hamilton-Jacobi equation is 

dS ^ / dS\ ^ 

with S{q,qo,t). A canonical transformation (ps ■ (901P0) ~^ il^P) is 
generated using 

dS dS 
dqo oq 

The initial conditions are chosen so that S generates an identity trans- 
formation at t = 

The generating function is expanded 

°° At" 
S = ^2 — 7-S.,{po,q) 

i=0 

and truncated to get an approximate solution Sn- The functions Si 
are computed recursively and they satisfy the reversibility condition for 
canonical Hamiltonian systems with the Hamiltonian H = ^p^ + V(q) 

Si{q, qo.,t) = ~Si{qo, q, -t) 

The Ge-Marsden algorithm in 24, is based on a reduction technique 
which is in fact a symplectic integrator on the symplectic leaves. The 
resulting scheme preserves the Lie-Poisson structure on the reduced 
space, i.e. it is symplectic on the symplectic leaves and preserves the 
angular momentum. Ge-Marsden algorithm is an implicit scheme for 
general Lie-Poisson systems, but it becomes explicit for the Euler top 
equations The application of the Ge-Marsden algorithm for the 
Euler top equations lO was also given in [23 ■ Similarly in gen- 
erating functions are developed for general Poisson maps, which are 
symplectic on the symplectic leaves. The Ge-Marsden algorithm is 
then developed using exponential maps in the general Lie algebra set- 
ting in and applied to truncated finite dimensional Vasto-Poisson 
equations. Numerical results show the long term preservation of the 
energy and therefore the stability of the algorithm. The Ge-Marsden 
algorithm is composed to obtain higher order Lie-Poisson integrators 
in |7| and applied again to Vasto-Poisson equations. A modification 
of the Ge-Marsden algorithm is given in 31 and this is illustrated for 
the Euler top equation 0. 
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5.2 Lie-Poisson integrators 

The implicit midpoint scheme preserves the constant Poisson structure 
[15] . It was shown in |S| that the imphcit midpoint scheme is an almost 
Poisson integrator, that it preserves the Lie-Poisson structure up to the 
second order terms of time step size At. The invariance of symplec- 
tic difference schemes under symplectic transformations is discussed 
in [22]. The Hamiltonian function takes different forms in different 
coordinates, so the symplectic schemes are coordinate dependent, i.e.. 
they are not covariant. This can be shown for the midpoint scheme, 
which is not invariant under general symplectic coordinate transforma- 
tions. Let the midpoint scheme for the canonical Hamiltonian system 
ii = J^^Hu be 

Un+l - Un 1 fl, \ , . 

— = J Hu ( 2 + ) 

In an another coordinate system v, u — S{v), where S' is a symplectic 
coordinate transformation, the Hamiltonian becomes H(S{v))eLiid the 
midpoint scheme becomes 

" - J-'H.. ( Sihv„ + v„+,)] (20) 



At " V 2 

while the original midpoint scheme H19|) in the coordinates v is 

S{Vn+l) - SjVn) 1 A/o/ ^.Q^ \\ rn^\ 

= J l-{S{v„) + S{Vn+l) j (21) 

If S is nonlinear then the schemes H2UI) and (|21|l are not identical. If 
5 is a linear symplectic transformation JS — J, then (|2U|) and (|21|l 
are identical. It was shown that some symplectic difference schemes 
are invariant under certain symplectic coordinate transformations and 
such schemes are studied for Lie-Poisson systems in |22| . 

For Lie-Poisson systems in three dimensions, it is more convenient 
to work with the one form j25| . For the Euler top equation the one 
form to be preserved becomes ^S] 

3 

J = ^ = eijkJtjdxk = 2'^Xidxi (22) 

It was shown in ^H) that 1)22(1 is preserved by the symplectic Runge- 
Kutta methods. Because the energy integral (7) and the total angular 
momentum (8) of Euler top equation are quadratic, they are preserved 
exactly by symplectic Runge-Kutta methods of Gauss-Legendry family 
and by partitioned Runge-Kutta methods of Lobated ILIA-B type. A 
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detailed comparison of various methods for the rigid body equations 
can be found in [TT) . 

ExpHcit Lie-Poisson integrators can be constructed, if the Hamil- 
tonian is separable as for the canonical Hamiltonian systems. The 
linear structure of flows of Lie-Poisson systems were exploited for the 
construction of explicit Lie-Poisson integrators and in 0^ for 
the Euler top equation, truncated Vasto-Poisson equations, and Eu- 
ler equations in fluid dynamics. If the Hamiltonian can be written as 
H = X]"=i Hi{xi) , then each vector field Xi = JVHi can be integrated 
separately by an integrator (j)i{At). In the case of the Euler top equa- 

tion (joj one obtains, with the Hamiltonians Hi = 27-, i — 1, . . . , 3, 
the following splitting: 

ii = 

• X1X3 

• _ X1X2 

X3 — 

A composition of the discrete flows 

M^t) ° ■ ■ ■ ° M^t) 

gives a first order explicit integrator for the whole system. A reversible 
second order composition can also be constructed 

(/)i(At/2) o • • • o (j)niAt) o • • • o (j)i{At/2). 

The Casimirs, which are the angular momenta of the Euler top equa- 
tion, are preserved exactly. The Hamiltonian splitting technique was 
applied to some Lie-Poisson systems, which are obtained as finite di- 
mensional models of some Hamiltonian pde's in plasma physics and 
fluid dynamics. The truncated Vasto-Poisson equations and the sine- 
Euler equations [59[ as finite dimensional approximations to the two- 
dimensional Euler equation were integrated in |34j using the Hamilto- 
nian splitting technique and it was shown that all Casimirs are pre- 
served up to the round-off errors. It turns out that explicit Lie-Poisson 
integrators are by far faster than the Lie-Poisson integrators using gen- 
erating functions in |12l . 

Time reversible integrators of arbitrary order were constructed in 
for some three dimensional bi-Hamiltonian systems like Lotka- 
Volterra equations, Lorenz equation and for the Toda lattice. Be- 
cause the Hamiltonians of all these systems are separable for both 
Poisson brackets, methods based on Hamiltonian splitting can be ap- 
plied. Symplectic composition of split flows gives an explicit second 
order time-reversible integrator. Higher order composition techniques 
are used to increase the order of the methods in [23 ■ The Hamil- 
tonians are preserved over long time intervals and the periodicity of 
solutions were retained. 



Xi = 
X2 = 
±3 = 



_ X2X3 


X-LX2 
I2 



Xi 
X2 

is 



X2X3 

h 

_ X1X3 

h 
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In |20| for nondependence Poisson systems with 2n-dimensional 
symplectic structure K{x) — J~^{x) symplectic difference schemes 
were constructed using a parameterized nonUnear transformation from 
the Poisson manifold to the canonical symplectic space. But as stated 
in |20| . in general it is difficult to find parameters which preserve the 
Poisson structure. 

In l^nj, the integration of Poisson systems with constant Poisson 
structure is considered and it was concluded that from the class of 
symplectic Runge-Kutta methods, only the diagonally implicit ones 
preserve the constant Poisson structure. However Mclachlan pointed 
out in that because all Runge-Kutta methods are equivalent under 
linear change of coordinates, the constant structure is also equivalent 
to a canonical Poisson structure, therefore all symplectic Runge-Kutta 
methods preserve the constant Poisson structure. 

5.3 Semidiscretized Hamiltonian pde's and integrable 
discretizations 

Nonlinear pde's in the Hamiltonian form can be discretized in space, 
so that the resulting ode inherits the integrals of the pde, this process 
is known as the integrable discretization. There are some examples of 
nonlinear pde's like the Landau-Lifschitz equation and the NLS equa- 
tion, which possess integrable discretized nonlinear lattice equations. 

Example 8: Landau Lifschitz equation 

The Landau-Lifschitz (LL) equation is a Lic-Poisson Hamiltonian pde 

[HEI] 

— = S xV'^S + S X DS (23) 

and was proposed originally as a model of an anisotropic Heisenberg 
ferromagnet, where S = S{x,y,t) £ M^, ||S'||2 — 1 denote the spin 
length, the matrix D G M.^'^^ represents the anisotropy and may be 
assumed diagonal with |c?i| < \d2\ < Ida]. 

The semi-discretized form of LL on a one-dimensional domain = 
{—L, L) gives with {xi = —L + ih, i — 0, . . . ,n} the one-dimensional 
isotropic Heisenberg spin chain 

S^ = ^S^x{S,^, + S,+i). (24) 

The two-dimensional discretized equation over a rectangular domain 
fi = {—L,L) X {—L,L) with {{xi,yj) = {—L + ih,—L + jh),i,j = 
0, . . . , where h = 2L/n and n a positive even integer, has the form 

Sij — Sij X AI{Sij-i + Sij+i + Si-i,j + Si+i,j) (25) 
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where 

Sij « S{ih,jk,t) — — (Sij-i + 5'ij+i + Si-ij + Si+ij) 

and M = I/h'^ + D/A. As boundary conditions periodic or homogenous 
boundary conditions are taken. The equations (|24|) and H25|l are Lie- 
Poisson systems and are referred to, for h — 1, as lattice Landau 
Lifschitz (LLL) equations. In two dimensional case, the vector S is 
defined by using the natural ordering on the grid points: 

S — [Sf i, . . . , sf j^ sJi, . . . , • • • S'^i, . . . , S'Jn]'^ 

the equation (25) can be written as a Lic-Poisson system 

S = J{S)\7sH 

with the Hamiltonian 

H = — SijAI{Sij-i + Sij+i + Si-ij + Si+ij). 

id 

We use the following notation to describe the form of the structure ma- 
trix J{S). For a vector u e R'^, we associate a 3 x 3 skew-symmetric 
matrix skew(u) such that, for any other vector v G R'^, v x u = 
skew(u)?;. Then the structure matrix J{S) consists of blocks of the 
skew-symmetric matrices skew(S'ij) on the diagonal. The LL equa- 
tions (24) and (25) have two additional integrals; the individual spin 
lengths 115*^112 and the total spin length J^i ^i- This system was inte- 
grated for the one dimensional case in |23 by the splitting of the total 
Hamiltonian H = — Si ■ >5'.i+i a,s H = Hi + H2 

i?i - ^ 5i • 5*4+1, H2 = - ^^Si ■ Si-i, for i odd 

i i 

The resulting Lie-Poisson systems 

5' = J(5)Vffi, S=J{S)VH2 

are then integrated with a Lie-Poisson integrator and the corresponding 
flows are composed in a symmetric form. Another discretization is 
based on the partitioning of the whole lattice equation. In case of the 
one dimensional semidiscrete system the spins are ordered as 

Pi = S2i~i, Qi = 821, i = l,...,n/2. 

We obtain then the following equations 

P^y< {Q^-l+Q^), Q^ ^ Q^ X {P^ + P^+l) , i = 1, . . . , 7l/2. 
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Defining a vector filed splitting Vh — Vi + V2 by 



Qi X (Pi + P2) 


Q, X (P, + P,+i) 


Qn/2 X (P„/2+Pl) . 

one can easily see that each of the vector fields are exactly integrable. 
This system was integrated in |21| by a semi-explicit staggered scheme 
which correspond to the second order Lobatto IIIA-B pair |2()j . 

Both methods are time-reversible and preserve the individual lengths 
of the spin vectors and total spin length. In terms of operation counts 
the staggered scheme is cheaper than the Lie-Poisson integrator, espe- 
cially in two and three dimensions. 



Pi X {Qn/2 + Ql) 





Pi X (Q,_i + Q,) 




Vo = 



Pn/2 X (Q(„_l)/2 + Qn/2) 





Example 9: Nonlinear Schrodinger equation 

The semi-discretization of the NLS equation Q by finite differences or 
by spectral methods and the time integration of the resulting integrable 
and nonintegrable ode's was the subject of many papers in recent years. 
The standard finite discretization scheme results in the direct NLS 
system (DNLS) J 

= - 2C/, + C/,+i) - 2\Uj\^Uj (26) 

where h = L/N denote the grid spacing. DNLS (|26|l possesses canon- 
ical Hamiltonian structure with the Hamiltonian 




and with the only one additional integral of motion (the energy) / = 
J2j=o^ l^jP- This system has only two integrals whereas the NLS has 
infinitely many integrals, therefore the DNLS (|26|) is not integrable. An 
integrable discretization of Q is given by the Ablowitz Ladik discrete 
NLS equations (IDNLS): 

lil, = -^{Uj-i - 2U, + C/,+i) - |C/,f (C/,-1 + C/,+i) (27) 
with the Poisson bracket 

{qm,Pn} = 7-(l + qnPn)5m,n, {qrn,qn} = {Pm,Pn} 
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the Hamiltonian 



. N-l 

H = ~^Y1 + gn+l) + 2 log(l + h^qnPn)) 

3=0 

and an additional first integral 

N-l 

^=Y1 1i(Pi-^ +Pj+i) 
j=o 

where Pn — Un and = U*. The equation H27|l has n integrals, 
i.e. it is completely integrable. The nonintegrable DNLS (|26|l was 
integrated in (^7] by symplectic Runge-Kutta methods, which preserve 
the Hamiltonian and the energy integral. The integrable discretization 
of NLS, the IDNLS equations (|27|l were transformed in the canonical 
symplectic form using a Darboux transformation in 57 . The resulting 
system was then integrated with the symplectic Rungc-Kutta methods 
which preserve the Hamiltonian and the integrals in long-term, whereas 
the nonsymplectic methods fail. 

The integrable discrete NLS equation H27|) was integrated in 
by the construction of a generating function to preserve the Poisson 
structure. A linear drift in the Hamiltonian error was observed, and it 
was eliminated by a modification of the scheme in j2Z| ■ 

The following splitting of the Ablowitz-Ladik discrete NLS equa- 
tions was introduced in |54| : 

F--flow i^"-fiow F+-fiow 

q'n^ qn-li^-qnPn) -2(7„ +9,1+1 (1 - (ZriPn) 

P>7i= -Pn+l{l - qnPn) +2pn -p„- 1 ( 1 - (7„p„) 

The flows and are integrated by the symplectic Euler method 
and F^ by the implicit midpoint rule and the whole scheme is ob- 
tained by a symmetric decomposition. It was also shown that the 
Poisson structure for each flow was preserved by the symplectic Euler 
method O] . Unfortunately the preservation of the Poisson structure 
for numerical methods is not independent of the particular equations 
as in case of the Hamiltonian systems with symplectic structure. A 
method which is Poisson preserving for one system may not preserve 
the Poisson structure of another system |^. The integration above 
requires, unfortunately, very small step sizes which makes it useless 
for numerical computations. 

Example 10: Volterra lattice equations: 

The Volterra lattice equations (|16|l can also be represented with Ui = 
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j/2i-i, Vi — y2i in partitioned form 

Ui = Ui{v^ - v^-i), Vi = Vi{ui+i - Ui), i = l,...,m/2 (28) 

for even m [SS]. The equations (|28|) are also bi-Hamiltonian with the 
quadratic and cubic Poisson brackets hke the original Volterra lattice 
equations. They were integrated with the partitioned Lobatto IIIA- 
B methods in (30]. It was shown that the Poisson structure is pre- 
served by the first order symplectic Euler method and by the second 
order Lobatto IIIA-B methods, whereas any of the Poisson brackets 
of the Volterra lattice IjKil) can be preserved by the implicit mid-point 
method. 

Several geometric integrators for the NLS were compared recently 
in among them the multisymplcctic integrators. An alterna- 

tive formulation of the Hamiltonian pde's involves a local concept of 
symplecticness in space and in time . Many Hamiltonian pde's are 
integrated using finite difference and spectral methods by multisym- 
plcctic integrators, which show better conservation properties than the 
symplectic integration of semi-discretized pde's (see for example jSOl)- 
An alternative formulation of multisymplecticness and the variational 
integrators, which are designed to integrate the Hamiltonian pde's are 
given in • These methods also have good preservation properties. 
Which geometric integrators are the best for Hamiltonian pde's, is still 
an open question. 

6 Nambu-Hamilton systems and discrete 
gradient methods 

In 1973 Nambu introduced a generalization of the classical Hamil- 
tonian mechanics and an extension of the Poisson bracket. In this 
formulation, the Poisson bracket is replaced by a ternary or n-ary 
operation, called Nambu bracket. The underlying idea of this new 
formulation was that in statistical mechanics the basic result is the Li- 
ouville theorem, which follows from but does not require the Hamilton 
dynamics. 

Definition 6 In the Euclidean space with the coordinates xi and 
X2 the Nambu bracket can be written as 

SF^ F2} ^ _ mm ^ d{F^,F2) 

dxi 8x2 3x2 dxi d{xi,X2) 

where this bracket satisfies the Jacobi identity The canonical 

Nambu bracket is defined for a triple of variables Fi, F2, Ft, in the phase 
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space M'^ with the coordinates xi,X2,X3 by the formula 

d[Xi,X2,X3) 

where the right-hand side stands for the Jacobian. 

The first Nambu formulation was given for the Euler top equation 
with the two Hamiltonians Hi, the energy integral (0 and H2, the 
total angular momentum (jH)): 

dx 

— ^{Hi,H2,x} = VHi XVH2. 
dt 

The Liouville theorem is still valid; the corresponding phase flow pre- 
serves the phase space volume. The orbits of a Nambu-Hamilton sys- 
tem in the phase space are determined as the intersection of the sur- 
faces Hi =const and H2 =const. For the Euler top equation, it is given 
by the intersection of the sphere described by the total angular mo- 
mentum H2, and the total energy Hi. The trajectories can be found 
by Jacobi elliptic functions. Nambu introduced also a n— dimensional 
generalization of the Euler top equation : 

dx£ _ n-^ dHi dHk dHn-i 

dt ^ ''^ dxj dxk dxk 

jk---l 

dx_ _ d{x,Hi, ■ ■ ■ 

dt d{xi, - ■ ■ ,Xn) 

where eijk - i denotes the Levi-Civita tensor. In 1994 Takthtajan 
made a geometric formulation of Nambu-Poisson bracket, which sat- 
isfies the skew symmetry, Leibniz rule and the fundamental identity, 
which is a generalization of the Jacobi identity. According to . 
the dynamics on the Nambu-Poisson manifold is determined by n — 1 
Hamiltonians Hi , . . . , i?n-i and is described by the generalized Nambu- 
Poisson equation: 

dx 

— = {Hi, . . .,Hn-l,x} 

The fundamental identity is given the n = 3 as 

{{G,H,Fi},F2,F3} + {Fi,{G,H,F2},F3} + {Fi,F2,{G,H,Fi,F3}} 
= {G,H,{Fi,F2,F3}. 

A function F is called an integral of motion for the Nambu-Hamilton 
system if its Nambu bracket with the Hamiltonians Hi, . . . , Hn-i van- 
ishes. A Nambu bracket of n integrals of motion is again an integral 
of motion, this shows that many integrable systems possess Nambu 
structure. 
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Example 11: Lagrange system: 

The so called Lagrange system [SU* which occurs in SU (2) monopoles 
on 

Xi=X2X^, ±2^X3X1, ±3^X1X2 

can be written in the Nambu form 

= {x,Hi,H2], Hi = ]^{x\ - xl), H2 = \^{x\ ~ Xg). 

It turned out that many physical systems can be written in Nambu 
form; like the Halpcn system |5()j . some vortex equations in fluid dy- 
namics in], some two dimensional incompressible flow equations |4;^j . 
SU{n) isotropic harmonic oscillators and the 5*0(4) Kepler problem 
. The three dimensional Lotka-Volterra equations (|12|l also can be 
written in Nambu form 

X = M{x)VHi X VH2 

where M(x) = X1X2X3 with the Hamiltonians Hi and H2 in (|13|l . 

Poisson systems can also be considered as a special case of a discrete 
gradient system 

X = S{x)VI{x) (29) 

where S{x) is skew-symmetric and I{x) is a first integral. The integral 
preserving methods based on the construction of discrete gradients 
in |37[ 117) are similar to those conserving schemes for Hamiltonian 
systems. The discrete gradient method consists of the replacement 
of the continuous gradient WI{x) by the discrete gradient VI{x){x,x) 
which satisfies 

yi{x){x, x){x -x) ^ V{x) - I{x), 

VI{x){x,x) = VI{x) +0{x ~ x). 

The resulting discrete gradient system can then be constructed as 

^^SVI{x)ix,x) 

where S is any consistent skew-symmetric matrix such that S{x, x) = 
S (^Y^)- It can be easily seen from the skew-symmetry of S, that the 
first integrals are preserved, i.e. I{x) = I{x). Discrete gradients are 
not unique; several examples of them are given in (371 1471 . 

Systems with multiple integrals Ii{x), ■ ■ ■ ,/„_i(x) like the Nambu 
systems can be written in the multi-gradient form 

X = S{x)VIi{x) ■ ■ ■ V/„_i(a;) 

where S(x) is an n skew-symmetric tensor. Discrete gradient systems 
which preserve the multiple gradients can be constructed similarly j^fij . 
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